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We perform a joint analysis of current data from cosmology and laboratory experiments to con¬ 
strain the neutrino mass parameters in the framework of bayesian statistics, also accounting for 
uncertainties in nuclear modeling, relevant for neutrinoless double P decay (Oz/2/3) searches. We find 
that a combination of current oscillation, cosmological and Oi/2/3 data constrains mpp < 0.045 eV 
(0.014 eV < mgp < 0.066 eV) at 95% C.L. for normal (inverted) hierarchy. This result is in practice 
dominated by the cosmological and oscillation data, so it is not affected by uncertainties related to 
the interpretation of Qu2P data, like nuclear modeling, or the exact particle physics mechanism un¬ 
derlying the process. We then perform forecasts for forthcoming and next-generation experiments, 
and find that in the case of normal hierarchy, given a total mass of 0.1 eV, and assuming a factor- 
of-two uncertainty in the modeling of the relevant nuclear matrix elements, it will be possible to 
measure the total mass itself, the effective Majorana mass and the effective electron mass with an 
accuracy (at 95% C.L.) of 0.05, 0.015, 0.02 eV respectively, as well as to be sensitive to one of the 
Majorana phases. This assumes that neutrinos are Majorana particles and that the mass mechanism 
gives the dominant contribution to 0v2P decay. We argue that more precise nuclear modeling will 
be crucial to improve these sensitivities. 

PACS numbers: 14.60.Pq, 23.40.-s, 98.80.-k 


I. INTRODUCTION 

It is by now firmly established by oscillation experiments 
that neutrinos do have a mass. However, oscillation ex¬ 
periments are only sensitive to neutrino mass differences 
and mixing angles, and thus do not provide information 
on the absolute scale of masses, on the mass hierarchy nor 
on their Dirac or Majorana nature. The nature of neu¬ 
trino masses and their smallness with respect to those 
of the charged leptons represents a puzzling fact, possi¬ 
bly related to the mechanism of neutrino mass genera¬ 
tion. Three main avenues are currently being pursued 
in order to experimentally probe the absolute scale of 
neutrino masses, namely i) direct measurements, study¬ 
ing the kinematics of /3 decay [T], ii) searches for neutri¬ 
noless double j3 decay (0:/2/3) [5], and iii) cosmological 
observations [S]. Approaches based on kinematic argu¬ 
ments have the advantage of being very direct and model- 
independent. An alternative is to study 0;/2/3 decay, i.e., 
the double /3 decay of nuclei, in which no neutrinos are 
present in the final state. If observed, it would guarantee 
that neutrinos have a non-vanishing Majorana mass [4j; 
if not, upper limits on the mass scale can still be placed, 
under the assumption that neutrinos are Majorana par¬ 
ticles. Relating the (potentially) observed rate for this 
process to neutrino masses also requires to assume that 
the mass mechanism is the dominant one leading to 0 j^ 2/? 
decay. It is worth noting that even if this is the most 
natural scenario, nevertheless other possibilities exist, in¬ 
volving additional physics beyond the standard model. 
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see e.g. Refs. [3 [6] for a discussion. Moreover, our 
imprecise knowledge of the appropriate nuclear matrix 
elements is a relevant source of uncertainty on the inter¬ 
pretation of the results of these experiments [2]- Finally, 
neutrino masses can be measured through cosmological 
observations, like measurements of the temperature and 
polarization anisotropies of the cosmic microwave back¬ 
ground, or of the distribution of large scale structures, 
since massive neutrinos affect the background evolution 
of the Universe, as well as the growth of cosmological 
perturbations. Cosmology presently provides the most 
stringent limits on the absolute scale of neutrino masses 
[7], with the shortcoming that these limits depend on 
assumptions on the underlying cosmological model. 

The three approaches outlined above should be seen 
as complementary, as each of them presents its own 
advantages and disadvantages, and also because they 
probe slightly different quantities related to the neutrino 
masses. For this reason, it appears natural to combine 
data from direct measurements, Qv2j3 searches and cos¬ 
mology, other than from oscillation experiments, in order 
to constrain the neutrino mass parameters [S]- In this 
paper, we want to derive joint constraints on neutrino 
mass parameters from the most recent observations from 
both laboratory and cosmological experiments, combin¬ 
ing them in the framework of Bayesian statistics. In 
particular, for Oz/2/3 experiments, we take into account 
the uncertainty related to nuclear matrix elements, by 
treating it as a nuisance parameter to be marginalized 
over, in order to account its impact on the neutrino mass 
estimates. We also perform forecasts, considering both 
forthcoming and next-generation experiments. 
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II. METHOD 


We use nii {i = 1, 2, 3) to denote the masses of the 
neutrino mass eigenstates Vi. We denote with 1 and 
2 the eigenstates that are closest in mass; moreover, 
we take m 2 > mi, so that is always positive, 

while the sign of discriminates between the nor¬ 

mal (NH) and inverted (IH) hierarchies, for Am§]^ > or 
< 0, respectively. The neutrino mass eigenstates are re¬ 
lated to the flavour eigenstates z/q, (a = e, fi, t) through 
Va = where Uai are the elements of the neu¬ 

trino mixing matrix U, parameterized by the three mix¬ 
ing angles (0i2, 623 , 613 ), one Dirac (<5) and two Majorana 
(<321, <331) CP-violating phases. Oscillation phenomena 
are insensitive to the two Majorana phases, that how¬ 
ever affect lepton number-violating processes like Oz^2/3 
decay. The different probes of the absolute scale of neu¬ 
trino masses are sensitive to different combinations of 
the mass eigenvalues and of the elements of the mix¬ 
ing matrix. /3 decay experiments measure the squared 
effective electron neutrino mass m'^ = 
while 0v2l3 searches are sensitive to the effective Ma¬ 
jorana mass m^/3 = where (j )2 = <321 and 

4>3 = 031 — 25. Finally, cosmological observations probe, 
at least in a first approximation, the sum of neutrino 
masses mi = mi + m 2 + m 3 . 

We perform a bayesian analysis based on a Markov 
chain Monte Carlo (MCMC) method, using cosmoMC 
[S] as a generic sampler in order to explore the pos¬ 
terior distribution of the parameters given the data. 
We consider the following vector of base parameters: 
[M^, Amh, Amh, sin^0i2, sin^dia, 4 >i, 4 ) 2 , (,) where ^ 
is a “nuisance” parameter related to the uncertainty in 
nuclear modeling (see below). We assume uniform prior 
distributions for all parameters. We do not consider the 
mixing angle 623 since none of the mass parameters de¬ 
pend on it. 

We consider data from oscillation experiments, direct 
measurements of the electron neutrino mass, 0 v 2 P 
searches and cosmological observations, all folded in the 
analysis through the corresponding likelihood function. 
Our baseline dataset is the most recent global fit of the 
neutrino oscillation parameters nni, updated after the 
Neutrino 2014 conference. We model the likelihood as a 
the product of individual gaussians in each of the oscil¬ 
lation parameters, since correlations can be neglected for 
our purposes [TUHIS- For the means and standard de¬ 
viations, we take respectively the best-fit value and the 
1(7 uncertainty quoted in Tab. II of Ref. uni. When 
the error is asymmetric, we conservatively take the stan¬ 
dard deviation equal to the largest between the left and 
right uncertainties. For direct measurements, we consider 
KATRIN [13] and HOLMES [T3] as our forthcoming and 
next-generation datasets, respectively. KATRIN is ex¬ 
pected to reach sub-eV sensitivity in mp, while HOLMES 
could go down to 100 meV. Kinematic measurements are 
directly sensitive to the square of the effective electron 
neutrino mass, so in both cases we take the likelihood to 


be a gaussian in (with the additional condition that 
m'p > 0), with a width given by the expected sensitivity 
of the experiment, i.e. cF{m'p) = 0.025, 0.006 eV^ for KA¬ 
TRIN and HOLMES, respectively. For 0v2P searches, we 
consider the current data from the GERDA experiment 
[Is] as the present dataset, its upgrade to the so-called 
“phase 2” for the near-future, and the nEXO experiment 
[IB] as a next-generation dataset. Oz/2/3 experiments are 
sensitive to the half-life of 0zz2/3 decay Assuming 

the Majorana nature of neutrinos, and that 0z^2/3 decay is 
induced by the exchange of light Majorana neutrinos (in 
the following we shall always assume that this is the case, 
unless otherwise stated), is related to the Majorana 
effective mass through: 


rpOl/ _ 

^ 1/2 — 


1 



( 1 ) 


where me is the electron mass, is a phase space fac¬ 
tor and is the nuclear matrix element. The phase 
I of the GERDA project provides the tightest bounds 
on the half-life of Oz/2/3 decay of ^®Ge, reporting a limit 
r0)^2 > 2-1 X 1025 yr at 90% G.L. {mpp < 200-600meV) 
BSQ The upgrade to the phase H of the experimental 
program is expected to increase the 90% G.L. sensitiv¬ 
ity to > 1-5 X 1025yr, {mpp < 90 — 150 meV) for 
40 kg of detector mass and 3 years of observations [19] . 
nEXO is a next-generation ton-scale experiment for the 
detection of 0 v 2 p decay of ^^^Xe, conceived as scaled-up 
version of the currently ongoing project EXO, with an 
estimated sensitivity > 6.6 x lO^’^yr at 90% G.L. 
{mpp <7—18 meV) for 5 tons of material and 5 years 
of data [30]. We model the likelihood of 0d2P exper¬ 
iments as a Poisson distribution in the number of ob¬ 
served events in the “region of interest” (the energy win¬ 
dow around the Q-value of the decay) with an expected 
value A = Ag -I- given by the sum of signal (S) and 
background (B) contributions. For a given value of T^j 2 ^ 
the expected number of signal events observed in a time 
Fobs for a detector mass M is 


= 


In 2 Nji£e 


TTlenr 


'T'Oi' ’ 
-^ 1/2 


( 2 ) 


where Na is Avogadro’s number, £ = MTohs is the 
exposure, e is the detector efficiency, menr is the mo¬ 
lar mass of the enriched element involved in the de¬ 
cay. The level of background is usually expressed in 
terms of the “background index”, i.e. the number of 
expected background events per unit mass and time 
within an energy bin of unit width. For GERDA-I, we 


^ Other isotopes currently yield T^j 2 > 4.0 X 10^"^ yr (^®'^Te) [l6| 
and ^ X 10^® yr (^®®Xe) at 90% C.L. ^Tji corresponding 
to mpp < 270 —760meV and mpp < 140 —280meV, respectively. 
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use the parameters reported in Tab. I of |15j for the 
case with pulse-shape discrimination. For GERDA-II, 
we consider a reduction of the background index down 
to 10“^ countskeV“^kg“^yr“^, a total exposure of 120 
kg yr, and the same efficiency as GERDA-I [^. For 
nEXO, we assume a background index corresponding to 
3.7 events ton“^yr“^ in the region of interest and an ex¬ 
posure of 25 ton yr [501 1 the same efficiency as EXO 
[55] . We also consider an update to nEXO in which the 
background in the inner 3 tons of the detector can be 
reduced by a factor 4 through Ba tagging. We assume 10 
years of observations for this updated version [50] . 

In order to account for the uncertainty related to nu¬ 
clear modeling |23j . including both that on nuclear ma¬ 
trix elements (NME) and that on the axial coupling con¬ 
stant, we compute for a given using fiducial 
values of these quantities, and then rescale it by a fac¬ 
tor ^. A similar approach was used in Ref. [51] in a 
frequentist framework, while we refer to Ref. [25] for a 
different bayesian approach. The fiducial values are gA = 
1.273 for the axial coupling, = 2.363 x 10“^^yr“^ 
(14.58 X lO-i^r"^) and MP'' = 3 (2) for ^^Qe (^^eXe). 
The value of ^ is extracted at every step of the MG from a 
uniform distribution in the range [0.5, 2], and marginal¬ 
ized over. This is equivalent, for example, to assume that, 
given exact knowledge of the axial coupling, the numeri¬ 
cal estimates of the NME can be wrong by up to a factor 
two in either direction. Finally, for what concerns the 
cosmological dataset, we use results obtained combining 
full mission Planck temperature and polarization data 
with data on the baryon acoustic oscillations |7] , as both 
our current and forthcoming reference dataset. For sim¬ 
plicity, we shall refer to this dataset simply as “Planck 
2015”. In particular, we use the chains publicly avail¬ 
able through the Planck Legacy Archive [26] to derive 
the posterior distribution of My given these data, cor¬ 
responding to a 95% upper limit My < O.I7eV. As a 
next-generation experiment, we consider the Euclid mis¬ 
sion. The combination of all Euclid probes (weak leas¬ 
ing tomography, galaxy clustering and ISW) with data 
from Planck is expected to constrain the sum of neutrino 
masses with a sensitivity of 0.06 eV for My = O.IeV, as 
reported in Tab.2 and the main text in ^7\. We shall 
refer to this dataset simply as “Euclid”. We model the 
likelihood as gaussian in My, with a {My) = 0.06 eV and 
the addition of the physical prior My > 0. 

To summarize, we consider four combinations of 
datasets. All of them include the most updated informa¬ 
tion from oscillation experiments. The “present” dataset 
includes Planck 2015 for cosmology, and GERDA-I for 
0u2l3 searches We do not include information from avail¬ 
able direct measurements (e.g. those from the Troisk 
and Mainz experiments) since they do not add infor¬ 
mation on 171/3 with respect to the data already consid¬ 
ered. The “forthcoming” dataset consists of the same 
cosmological data as the previous dataset, GERDA-II, 
and KATRIN for kinematic measurements. The “next 
generation I (II)” dataset includes Euclid, nEXO with¬ 


out (with) Ba tagging and HOLMES. For future data, 
we have to assume fiducial values of the parameters: in 
the case of the forthcoming dataset, we take them equal 
to their best estimates from the combination of oscilla¬ 
tions and Planck2015. For the futuristic case, we assume 
My = 0.1 eV and estimate rri/s and m/ 3/3 from the combi¬ 
nation of Euclid and oscillation parameters. 


III. RESULTS 

We present our results for My, mp and m/ 3/3 in Tab. [Jfor 
the three datasets described above. We report limits R)th 
in the case where ^ is fixed to 1, and when ^ is marginal¬ 
ized over, in order to show the impact of uncertainties in 
nuclear modeling. We quote our results, both in the text 
and table, in terms of the Bayesian 95% minimum credi¬ 
ble interval [55] . When this interval includes the minimal 
value of the parameter allowed by oscillation measure¬ 
ments, we only quote the extremes of the range; on the 
contrary, we report the mean ± the 95% uncertainty. We 
do this in order to emphasize a “detection” scenario - i.e., 
one in which the observations point to a value of the pa¬ 
rameter under consideration being different, with a given 
statistical significance, from the lowest value allowed by 
oscillations alone - from a “non-detection” scenario in 
which this oscillation minimal value is still allowed. We 
choose to identify the minimal value allowed by oscilla¬ 
tions as the bayesian 95% C.L. lower limit of the neu¬ 
trino mass parameters when the lightest eigenstate is set 
to zero. With this definition, we get M™'" = 0.057, eV 
(M™ = 0.096eV), = 0.009eV (m^“ = 0.047eV), 

= 0.002eV = 0.016eV) for NH (IH). We 

would like to point out that the exact definition of the 
minimal value is somehow arbitrary, in a sense that it is 
not formally well-defined due to the finite precision of the 
oscillation measurements. For example, we could have 
chosen the lowest value allowed by fixing the oscillation 
parameters to their best-fit values, rather than comput¬ 
ing the bayesian 95% lower limit, and this would equally 
make sense. This choice only affects the way in which 
limits are reported in Tab. [l] (and we verified that it also 
has a minor impact in that respect), so it does not alter 
our conclusions in any way. In any case, we recall that the 
confidence intervals represent a compression of the infor¬ 
mation contained in the one-dimensional posteriors, that 
fully represent the probability distribution associated to 
a given parameter. In Fig. [ijwe show the marginalized 
one-dimensional posterior distributions for the mass pa¬ 
rameters. In most cases, the low mass region is excluded 
by the oscillation data, with the only exception of m/ 3/3 in 
the case of NH; the reason is that in this case the phases 
can arrange in order to yield mpp = 0 even for finite 
values of the mass differences. Present data provide sim¬ 
ilar limits independently of whether nuclear uncertainties 
are marginalized over. This happens because the present 
constraints are dominated by the cosmological limit on 
My, that translates directly to bounds on mp and mpp 
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once oscillation data are taken into account (this can be 
understood by noticing that the direct limits on these pa¬ 
rameters are much weaker). We have verified explicitly 
that this is the case by performing parameter estimation 
using only Planck2015 and oscillation data, as shown in 
Figj^ In particular, we find that the present data con¬ 
strain m/ 3 ^ < 0.045 eV (0.014 eV < mpp < 0.066 eV) 
at 95% C.L. for NH (IH), regardless of the inclusion of 
Qv2P information. Forthcoming datasets yield similar 
constraints for the mass parameters; this means that the 
improved sensitivity of GERDA-II and the inclusion of 
KATRIN add only marginally to the Planck2015 plus os¬ 
cillations data combination. The fact that present and 
forthcoming limits on mpp are dominated by the lat¬ 
ter dataset has the consequence that they do not de¬ 
pend on the modeling of 0^2/3-decaying nuclei, nor on 
assumptions about the mechanism that induces the de¬ 
cay (while, on the other hand, they are affected by the 
model dependence of the cosmological analysis). This 
picture changes substantially for next-generation experi¬ 
ments. In this case, cosmological observations and {)v2fi 
searches have comparable power in constraining the mass 
parameters, and the nuclear uncertainties - as well as the¬ 
oretical assumptions about the particle physics of Qu2j3 
decay - play a role in deriving parameter constraints. We 
find that, if neutrinos are Majorana and 0^'2/3 decay is 
dominantly induced by the mass mechanism, marginal 
95% evidence for non-minimal mass parameters can be 
obtained in the case of normal hierarchy, even when nu¬ 
clear uncertainties are taken into account. This detection 
is further strengthened in the Next generation-II dataset, 
for which we get M^, = 116^51 meV, mpp = 20tj5meV, 
mp = 30 ± 20 meV. In the case of inverted hierarchy, 
we obtain upper limits for M^, and mp, and a more 
than 95% evidence for non-minimal mpp. In particu¬ 
lar, for the Next generation-II dataset and marginaliz¬ 
ing over nuclear uncertainties, we find < 187 meV, 
TTT-pP = 43i*i24 meV and mp < 69meV. Finally, we report 
that while present and forthcoming experiments have lit¬ 
tle, if none, sensitivity to the neutrino mixing phases, the 
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combination of next-generation experiments will possibly 
allow to determine the value of 021 , as shown in Fig. 


IV. CONCLUSIONS 

The combination of current and forthcoming data from 
oscillation, kinematic, Oi/2/3 and cosmological experi¬ 
ments allows to put upper bounds My < 0.19 (0.21) eV, 
mpp < 0.04 (0.06) eV and mp < 0.06 (0.08) eV for NH 
(IH). These limits are dominated by the combination of 
oscillations and cosmological data and as such are not af¬ 
fected by uncertainties in nuclear modeling, nor rely on 
the knowledge of the particle physics mechanism lead¬ 
ing to 0v2[3 decay. If neutrinos are Majorana particles 
and 0i'2(3 decay is induced by the exchange of light Ma¬ 
jorana neutrinos, and further assuming a total mass of 
0.1 eV and a factor 2 uncertainty in nuclear modeling, 
next-generation experiments will ideally allow to mea¬ 
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duced by up to a factor 4 by a better modeling of the 
nuclear factors. Next-generation experiments will also 
be sensitive to the phase 021 - 
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FIG. 1. Posterior distributions for the neutrino mass parameters, for NH (top row) and IH (bottom row). Solid (dashed) curves 
correspond to marginalization over nuclear uncertainties (fixed fiducial values for nuclear parameters). Note that the “Next- 
generation I” and “Next-generation H” forecasted results have been derived by assuming a fiducial value of Mi, = (0.10±0.06)eV. 



FIG. 2. Two-dimensional posterior distributions for the 
neutrino mass parameters in the mpp — mp plane, for NH 
(red) and IH (blue), from the combination of oscillation and 
“Planck 2015” datasets. Gontours correspond to 95% C.L.. 
The dashed lines are the 95% C.L. upper limits on mpp from 
GERDA phase 1 within the range ^ = [0.5 — 2] (vertical) and 
on from KATRIN (horizontal). 
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Parameter Present Forthcoming Next-generation Next-generation II 

Normal Inverted Normal Inverted Normal Inverted Normal Inverted 

M^,Cfree [50 - 194] [91-217] [50 - 188] [91-213] 109t^® [95 - 186] 116^5^ [95- 187] 

[50 - 195] [91-217] [50 - 188] [91-212] 105^17 [94- 176] llOtg [94- 177] 

mi3f3,^iTee <45 [14-66] <42 [14-64] [1-35] 43^24 20^11 4312° 

m/3/3, ^ = 1 <45 [14 - 66] < 42 [14 - 63] ISl^g 39 ± 9 16 ± 9 38 ± 6 

m/ 3 , Cfree [5-59] [45 - 78] [5-57] [45 - 76] [8-48] [47- 68] 30 ± 20 [47- 69] 

m/3,^ = l [5-59] [45 - 78] [5-57] [45 - 76] [8-45] [46- 65] 28^.1° [46 - 66] 

TABLE I. Limits on nentrino mass parameters from different datasets. For each parameter, we quote the 95% bayesian 
minimum credible interval. When this interval includes the minimal value for the parameters computed as the lower limit at 
95% C.L. allowed by oscillation data, we only quote the extremes of the range; otherwise, we quote the 95% range around the 
mean value. Units are in meV. Note that the “Next-generation” and “Next-generation H” forecasted results have been derived 
by assuming a fiducial value for the sum of the neutrino masses of Mi, = 0.10 eV. 



